Method for tracking the position and the heading of a vehicle using dead reckoning and a tracking device for carrying out the method

ABSTRACT

Position (x ref ) and heading (H) of a vehicle are tracked using wheel tick numbers coming from wheel tick units at the front wheels only from which a velocity (v) and a heading rate (h) are calculated. The calculations are based on an intermediate wheel angle (δ). An exact solution δ ex  for the same depends on the velocities of the left front wheel v l  and the right front wheel v r . An approximate solution δ app  is calculated according to a simpler formula and used wherever its deviation from the exact solution δ ex  is not significant, i.e., where the latter is not larger than a threshold wheel angle δ thr . Only for δ app &gt;δ thr  where the deviation is larger than a threshold, the wheel angle according to the exact solution δ ex  is used. The latter is determined by first establishing the approximate solution δ app  and reading the corresponding exact solution δ ex  from a look-up table.

RELATED APPLICATION

This application claims the benefit of priority under 35 USC §119 to EPPatent Application No. 10405212.1 filed on Nov. 4, 2010, the entirecontents of which are incorporated herein by reference.

FIELD OF THE INVENTION

The invention concerns a method for keeping track of the position andthe heading of a vehicle using tracking by dead reckoning. Such methodsare used in cars and other vehicles, often complementing determinationof the vehicle position by Global Navigation Satellite System (GNSS), inparticular, where GNSS signals are temporarily too weak or unavailable,e.g., in tunnels, car parkings and the like. The invention also concernsa tracking device for carrying out the method.

PRIOR ART

Methods for tracking the position and heading of a vehicle by deadreckoning using rear wheel velocities derived from sensor measurementsare well known. Although the evaluation of front wheel velocities ismore complicated due to the variability of the wheel angles methodsusing only front wheel velocities have also been proposed. So far onlyapproximate solutions connecting the velocity and heading rate of thevehicle to left and right front wheel velocities have been employed,see, e.g.: Ch. Hollenstein, E. Favey, C. Schmid, A. Somieski, D. Ammann:‘Performance of a Low-cost Real-time Navigation System usingSingle-frequency GNSS Measurements Combined with Wheel-tick Data’,Proceedings of ION GNSS 2008 Meeting, Savannah, September 2008, thecontents of which are incorporated herein in their entirety. Accordingto this paper, the wheel angle was determined by iterative methods. Theapproximations used can compromise the precision of the tracking wherelarge wheel angles cannot be excluded.

SUMMARY OF THE INVENTION

It is the object of the invention to provide a method for keeping trackof the position and heading of a vehicle which provides improvedprecision.

The advantages of the method according to the invention are particularlypronounced where large wheel angles occur, e.g., in underground carparkings and similar environments where, at the same time, GNSS signalstend to be weak or absent. The method can be implemented in ways whichare very economical as to processing power and memory requirements.

The sensor measurements can be used in combination with GNSSmeasurements where such are available, e.g., by feeding both to atightly coupled Kalman filter.

It is another object of the invention to provide a tracking devicesuitable for carrying out the method according to the invention.

BRIEF DESCRIPTION OF THE DRAWINGS

In the following the invention is described in more detail withreference to drawings which only illustrate an examplary embodiment.

FIG. 1 schematically shows a vehicle comprising a tracking deviceaccording to the invention and whose position and heading is beingtracked using the method according to the invention,

FIG. 2 shows a flow diagram of the method according to the invention,

FIG. 3 shows a part of the flow diagram of FIG. 2 in more detail, and

FIG. 4 shows a diagram with an approximate solution and an exactsolution for an intermediate wheel angle δ as shown in FIG. 1.

FIG. 5 shows a block diagram of one embodiment of the tracking device ofFIG. 1.

DESCRIPTION OF THE PREFERRED EMBODIMENTS

In vehicles so-called ‘dead reckoning’ is often used for tracking theposition of the vehicle by determining the velocities of the wheels andusing the results as a basis for updating the vehicle position, eitherin combination with GNSS or, where GNSS is not available, by itself.Dead reckoning and GNSS tracking calculations can be carried out in anappropriate tracking device, usually a component of the vehicle orattached to the same, and the results displayed or made use of in otherways.

A wheel velocity is normally determined from sensor measurements,usually signals known as ‘wheel ticks’ which are produced at particularrotational positions separated by fixed angular increments. Where sensormeasurements from the rear wheels are used, deriving the speed andheading rate of the vehicle is usually straightforward as wheel anglesmay be assumed to be constant and equal to zero.

In some cases, however, sensor measurements from rear wheels are notavailable and speed and heading rate are derived from front wheel sensormeasurements alone. This includes more involved calculations as thevariable wheel angles which are taken into account is derived from thewheel ticks as well unless they are available otherwise, e.g., fromdirect measurements.

In the following the case where wheel angles are derived solely fromfront wheel sensor measurements is described in detail. It is assumedthat the vehicle (FIG. 1) comprises a rigid frame 1 on which rear wheels2 a,b and front wheels 3 a,b are mounted. The vehicle is steerable inthat the front wheels 3 a,b are rotatable about vertical axes passingthrough their centers and may assume wheel angles covering an intervalcontaining 0° which value corresponds to a motion of the vehiclestraight ahead. A tracking installation comprises a tracking device 4which is fixed to the frame 1 and wheel tick units 5 a,b monitoring therotations of the front wheels 3 a;b or other sensors producing sensormeasurements from which the wheel velocities of the front wheels 3 a,bcan be derived. One embodiment of the tracking device 4 is shown in FIG.5. A GNSS antenna 6 which is also mounted on the frame 1, feeds, likethe wheel tick units 5 a,b, data to the tracking device 4. The trackingdevice 4 is configured to process the sensor measurements, i.e., signalsfrom the wheel tick units 5 a,b, in order to extract a position and aheading of the vehicle from them as well as GNSS signals received by theGNSS antenna 6 and contains the components for storing and processingdata and outputting results etc. which are used for the purpose. It isassumed that the frame 1 and the suspensions of the rear wheels 2 a,band the front wheels 3 a,b are symmetrical with respect to an axis ofsymmetry which passes through the midpoints between the rear wheels 2a,b and the front wheels 3 a,b.

In the following derivations it is assumed that there is no slippage andthat the momentary motion of the wheel center follows in each case thedirection implied by the wheel angle and, for simplicity, that thevehicle moves in a level plane defined by a northerly direction (n) andan easterly direction (e). Every vector x in the plane thereforeconsists of a northerly component and an easterly component, that is,x=(x^((n)), x^((e))). If the plane the vehicle moves on is inclined itis straightforward to project its motion onto a level plane.

In the example the state of the vehicle can be described by the positionof a reference point x _(ref)=(x_(ref) ^((n)), x_(ref) ^((e))) fixed onthe frame 1 of the vehicle and the heading H, i.e., the angle betweenthe northerly direction and the longitudinal direction of the saidframe, where the angle is measured clockwise.

In dead reckoning tracking of the vehicle position, a new state assumedby the vehicle at the end of a (k+1)^(th) time update interval of lengthΔt can be derived from the state at the end of the previous time updateinterval using equations of motion

$\begin{matrix}{x_{{ref},{k + 1}}^{(n)} = {x_{{ref},k}^{(n)} + {v\;{\cos\left( {H_{k} + \frac{h\;\Delta\; t}{2}} \right)}\Delta\; t}}} & (1) \\{x_{{ref},{k + 1}}^{(e)} = {x_{{ref},k}^{(e)} + {v\;{\sin\left( {H_{k} + \frac{h\;\Delta\; t}{2}} \right)}\Delta\; t}}} & (2) \\{H_{k + 1} = {H_{k} + {h\;\Delta\; t}}} & (3)\end{matrix}$where v is the momentary scalar velocity and h the momentary headingrate, i.e., the angular velocity of the vehicle. For every time updatestep the scalar velocity v and the heading rate h is thereforedetermined from the velocities v_(l) and v_(r) of the left and rightfront wheels. They correspond to raw velocities T_(l), T_(r) derivedfrom numbers of wheel ticks registered during the time interval Δtaccording tov _(l) =f _(l) T _(l)  (4a)v _(r) =f _(r) T _(r)  (4b)with wheel tick calibration factors f_(l) and f_(r) which depend on thewheel radii.

For simplicity it is (FIG. 1) assumed in the following that the vehicleis at the moment oriented in an easterly direction, i.e., that itsheading H equals 90°, and that the reference point x _(ref) of thevehicle, defined as the midpoint between the centers of the rear wheels,is at the origin O. Thereby the below derivations can be carried out ina simple manner in the fixed coordinate system described above. As theresults do not depend on the orientation of the vehicle and its positionon the plane this does not imply a loss of generality.

According to the above assumptions,x _(ref)=(0,0)  (5)and the midpoint between the centers of the front wheels is atx _(f)=(0,l)  (6)where l is the wheel base of the vehicle and the centers of the frontwheels are atx _(l)=(b,l)  (7a)x _(r)=(−b,l)  (7b)with b equalling one half of the front track gauge.

The vehicle, that is, its reference point x _(ref), moves at the momentwith a velocity v in an easterly direction:v =(0,v).  (8)

The corresponding velocities of the other relevant points x _(f), x_(l), x _(r) arev _(f)=(−v _(f) sin δ,v _(f) cos δ)  (9)andv _(l)=(−v _(l) sin δ_(l) ,v _(l) cos δ_(l))  (10a)v _(r)=(−v _(r) sin δ_(r) ,v _(r) cos δ_(r))  (10b)where v_(f), v_(l) and v_(r) each designate a scalar velocity, i.e., thelength of the respective vector and δ_(l), δ_(r) are the wheel angles ofthe left and the right front wheel, respectively, whereas anintermediate wheel angle δ is the angle between the easterly directionand the momentary direction of motion of x _(f).Defining

$\begin{matrix}{{a = \frac{b}{l}},{with}} & (11) \\{{\tan\;\delta}\; = \frac{l}{r}} & (12)\end{matrix}$where r is the distance between x _(ref) and the momentary center ofrotation of the vehicle C at (−r,0)and

${\tan\;\delta_{l}} = \frac{l}{r + b}$one gets

$\frac{1}{\tan\;\delta_{l}} = \frac{1 + {a\;\tan\;\delta}}{\tan\;\delta}$and finally

$\begin{matrix}{{\tan\;\delta_{l}} = {\frac{\tan\;\delta}{1 + {a\;\tan\;\delta}}.}} & \left( {13a} \right)\end{matrix}$Similarly, from

${\tan\;\delta_{r}} = \frac{l}{r - b}$one has

$\begin{matrix}{{\tan\;\delta_{r}} = {\frac{\tan\;\delta}{1 - {a\;\tan\;\delta}}.}} & \left( {13b} \right)\end{matrix}$From the rigidity of the vehicle it follows that

$\begin{matrix}{{\frac{\mathbb{d}}{\mathbb{d}t}{{{\underset{\_}{x}}_{f} - {\underset{\_}{x}}_{ref}}}^{2}} = {{2\left\langle {{{\underset{\_}{x}}_{f} - {\underset{\_}{x}}_{ref}},{{\underset{\_}{v}}_{f} - \underset{\_}{v}}} \right\rangle} = 0}} & (14)\end{matrix}$and together with (5), (6), (8) and (9) that l(v_(f) cos δ−v)=0 and

$\begin{matrix}{v_{f} = \frac{v}{\cos\;\delta}} & (15)\end{matrix}$and finally with (9)v _(f)=(−v tan δ,v).  (16)Using the same method for x _(l), x _(r) yields

$\begin{matrix}{v_{l} = \frac{v}{\cos\;{\delta_{l}\left( {1 - {a\;\tan\;\delta_{l}}} \right)}}} & \left( {17a} \right) \\{{v_{r} = \frac{v}{\cos\;{\delta_{r}\left( {1 + {a\;\tan\;\delta_{r}}} \right)}}}{and}} & \left( {17b} \right) \\{{\underset{\_}{v}}_{l} = {\frac{v}{1 - {a\;\tan\;\delta_{l}}}\left( {{{- \tan}\;\delta_{l}},1} \right)}} & \left( {18a} \right) \\{{\underset{\_}{v}}_{r} = {\frac{v}{1 + {a\;\tan\;\delta_{r}}}\left( {{{- \tan}\;\delta_{r}},1} \right)}} & \left( {18b} \right)\end{matrix}$and with (13a) and (13b), respectively, because of

$\begin{matrix}{\frac{1}{1 - {a\;\tan\;\delta_{l}}} = \frac{1}{1 - \frac{a\;\tan\;\delta}{1 + {a\;\tan\;\delta}}}} \\{= {1 + {a\;\tan\;\delta\mspace{14mu}{and}\frac{1}{1 + {a\;\tan\;\delta_{r}}}}}} \\{= \frac{1}{1 + \frac{a\;\tan\;\delta}{1 - {a\;\tan\;\delta}}}} \\{= {1 - {a\;\tan\;\delta}}}\end{matrix}$v _(l) =v(−tan δ,1+a tan δ)  (19a)v _(r) =v(−tan δ,1−a tan δ)  (19b)andv _(l) =v√{square root over (1+2a tan δ+(1+a ²)tan²δ)}  (20a)v _(r) =v√{square root over (1−2a tan δ+(1+a ²)tan²δ)}.  (20b)The scalar heading rate is

$\begin{matrix}{{h} = {\frac{\left( {v_{f} - v} \right)_{\bot}}{l} = \frac{{v_{f}^{(n)} - v^{(n)}}}{l}}} & (21)\end{matrix}$where the subscript ⊥ signifies the component of the vector which isperpendicular to the longitudinal axis of the vehicle, which leads to

$\begin{matrix}{h = \frac{v\;\tan\;\delta}{l}} & (22)\end{matrix}$where it is assumed that a heading rate in the clockwise directioncarries the plus sign.Defining V₊ and V⁻ asV ₊ =v _(l) ² +v _(r) ²  (23)andV ⁻ =v _(l) ² −v _(r) ²  (24)one gets from (20a) and (20b)V ₊=2v ²[1+(1+a ²)tan²δ]  (25)andV ⁻=4av ² tan δ  (26)or

$\begin{matrix}{{\tan\;\delta} = \frac{V_{-}}{4\;{av}^{2}}} & (27)\end{matrix}$and

$V_{+} = {{2{v^{2}\left\lbrack {1 + {\left( {1 + a^{2}} \right)\frac{V_{-}^{2}}{16a^{2}v^{4}}}} \right\rbrack}} = {{2v^{2}} + {\left( {1 + a^{2}} \right)\frac{V_{-}^{2}}{8a^{2}v^{2}}}}}$and finally

$\begin{matrix}{{v^{4} - {\frac{V_{+}}{2}v^{2}} + \frac{\left( {1 + a^{2}} \right)V_{-}^{2}}{16\; a^{2}}} = 0.} & (28)\end{matrix}$This quadratic equation for v² has a solution

$\begin{matrix}{v^{2} = {\frac{1}{4}\left( {V_{+} \pm \sqrt{V_{+}^{2} - {\frac{1 + a^{2}}{a^{2}}V_{-}^{2}}}} \right)}} & (29)\end{matrix}$with two branches, one with a plus-sign and one with a minus sign infront of the square root.With

$\begin{matrix}{P = \frac{V_{-}}{V_{+}}} & (30)\end{matrix}$(29) becomes

$\begin{matrix}{v^{2} = {\frac{V_{+}}{4}\left( {1 \pm \sqrt{1 - P^{2} - \left( \frac{P}{a} \right)^{2}}} \right)}} & (31)\end{matrix}$and for v one gets accordingly

$\begin{matrix}{{v = {{\pm \frac{1}{2}}\sqrt{V_{+}\left( {1 \pm \sqrt{1 - P^{2} - \left( \frac{P}{a} \right)^{2}}} \right)}}}\;,} & (32)\end{matrix}$where positive v signifies forward motion and negative v backward motionof the vehicle. (27) and (31) lead, together with (30), to exactsolutions of equations (20a), (20b) for the tangent of the wheel angle δand the said wheel angle itself:

$\begin{matrix}{{{\tan\;\delta_{ex}} = \frac{\frac{P}{a}}{1 \pm \sqrt{1 - P^{2} - \left( \frac{P}{a} \right)^{2}}}},} & (33) \\{\;{\delta_{ex} = {{\arctan\left( \frac{\frac{P}{a}}{1 \pm \sqrt{1 - P^{2} - \left( \frac{P}{a} \right)^{2}}} \right)}.}}} & (34)\end{matrix}$(33) can be inverted and P represented as a function of y=tan δ_(ex):

$\begin{matrix}{{P(y)} = {\frac{2{ay}}{1 + {\left( {1 + a^{2}} \right)y^{2}}}.}} & (35)\end{matrix}$This antisymmetric function has extrema at

$\begin{matrix}{{y_{\max} = \frac{1}{\sqrt{1 + a^{2}}}},} & \left( {36a} \right) \\{{y_{\min} = {- \frac{1}{\sqrt{1 + a^{2}}}}}{where}} & \left( {36b} \right) \\{{P_{\max} = {{P\left( y_{\max} \right)} = \frac{a}{\sqrt{1 + a^{2}}}}},} & \left( {37a} \right) \\{P_{\min} = {{P\left( y_{\min} \right)} = {- {\frac{a}{\sqrt{1 + a^{2}}}.}}}} & \left( {37b} \right)\end{matrix}$These are the values of P where the argument of the square root in (31)is zero and δ_(ex) assumes the values

$\begin{matrix}{{{\delta_{ex}\left( P_{\max} \right)} = {\arctan\left( \frac{1}{\sqrt{1 + a^{2}}} \right)}},} & \left( {38a} \right) \\{{\delta_{ex}\left( P_{\min} \right)} = {- {{\arctan\left( \frac{1}{\sqrt{1 + a^{2}}} \right)}.}}} & \left( {38b} \right)\end{matrix}$The branch of the solution (31) which corresponds to a range for δ_(ex)containing 0—motion of the vehicle straight ahead—is the first branchwhere the square root in (31) and consequently in (33), (34) is precededby a plus sign and which applies for 0≦|δ_(ex)|≦δ_(lim) withδ_(lim)=δ_(ex)(P_(max)) whereas the second branch, with a minus sign infront of the square root, covers δ_(lim)≦|δ_(ex)|≦90°. This implies thatif |δ_(ex)| may assume values larger than δ_(lim) the solution is noteverywhere unique as there are subranges of δ_(ex) surrounding δ_(lim)and −δ_(lim) where two solutions exist and it is not possible to pickthe correct one on the basis of the momentary wheel velocities alone.However, for δ_(ex)=δ_(lim) we have, by (13b),

$\begin{matrix}{\delta_{r} = {\delta_{wlim} = {\arctan\left( \frac{1}{\sqrt{1 + a^{2}} - a} \right)}}} & (39)\end{matrix}$and as long as the wheel angle of the right front wheel δ_(r) is, e.g.,by mechanical constraints, limited by δ_(r)≦δ_(wlim)—and, bysymmetry—the wheel angle of the left front wheel δ_(l) byδ_(l)≧−δ_(wlim)—the solution is unique as only the first branch in (31)and (33), (34) can apply. As δ_(wlim) is always larger than 45° and withmost vehicles the wheel angles are limited to values smaller than thatthe limitation 0≦|δ_(ex)|≦δ_(lim) will be assumed to hold in thefollowing.With this assumption (31) to (34) can be restated as

$\begin{matrix}{{v^{2} = {\frac{V_{+}}{4}\left( {1 + \sqrt{1 - P^{2} - \left( \frac{P}{a} \right)^{2}}} \right)}},} & (40) \\{{v = {{{\pm 1}/2}\sqrt{V_{+}\left( {1 + \sqrt{1 - P^{2} - \left( \frac{P}{a} \right)^{2}}} \right)}}},} & (41) \\{{\tan\;\delta_{ex}} = \frac{\frac{P}{a}}{1 + \sqrt{1 - P^{2} - \left( \frac{P}{a} \right)^{2}}}} & (42) \\{and} & \; \\{\delta_{ex} = {{\arctan\left( \frac{\frac{P}{a}}{1 + \sqrt{1 - P^{2} - \left( \frac{P}{a} \right)^{2}}} \right)}.}} & (43)\end{matrix}$The heading rate h can then be derived from (22), (41) and (42) as

$\begin{matrix}{h = {{\pm \frac{\sqrt{V_{+}}}{2\; b}}{\frac{P}{\sqrt{1 + \sqrt{1 - P^{2} - \left( \frac{P}{a} \right)^{2}}}}.}}} & (44)\end{matrix}$

The exact solutions derived above for v, δ and h are rather complicatedand, in particular, lead to unwieldy expressions for the derivativeswhich are needed where the dependence of v and h upon the scalar wheelvelocities v_(l) and v_(r) are used in a Kalman filter or a leastsquares or similar evaluation. Working with them in real time situationswhere processing power and time are limited is therefore often notpracticable. There are, however, more tractable approximations which areusually good enough, at least as long as the wheel angles are not toolarge. For large wheel angles the said approximations can then beappropriately corrected.

A relatively simple approximate solution which is well known per se fromgeometrical considerations can be developed starting with (20a) and(20b). (20a) can be rewritten as

$\begin{matrix}{v_{l} = {v{\sqrt{\frac{1}{\cos^{2}\delta} + {2\; a\;\tan\;\delta} + {a^{2}\tan^{2}\delta}}.}}} & (45)\end{matrix}$The square root can then be expanded, i.e., linearised, about

$\frac{1}{\cos^{2}\delta}$which yields

$\begin{matrix}{v_{l} = {\frac{v}{\cos\;\delta} + {v\;\cos\;\delta\; a\;\tan\;\delta} + {{O\left( \left( {a\;\tan\;\delta} \right)^{2} \right)}.}}} & (46)\end{matrix}$Disregarding the higher order terms, one arrives at

$\begin{matrix}{v_{l} = {\frac{v_{app}}{\cos\;\delta_{app}} + {{av}_{app}\sin\;\delta_{app}}}} & \left( {47\; a} \right) \\{v_{r} = {\frac{v_{app}}{\cos\;\delta_{app}} - {{av}_{app}\sin\;\delta_{app}}}} & \left( {47b} \right)\end{matrix}$where the suffix app signifies ‘approximate’ and v_(r) has been derivedin a manner entirely analogous to the derivation of v_(l) above.From (47a), (47b) follow

$\begin{matrix}{{v_{l} + v_{r}} = \frac{2\; v_{app}}{\cos\;\delta_{app}}} & (48)\end{matrix}$andv _(l) −v _(r)=2av _(app) sin δ_(app).  (49)Defining

$\begin{matrix}{p = \frac{v_{l} - v_{r}}{v_{l} + v_{r}}} & (50)\end{matrix}$one has

$\begin{matrix}{{p = {{a\mspace{11mu}\sin\;\delta_{app}\cos\;\delta_{app}} = {\frac{a}{2}\sin\; 2\;\delta_{app}}}}{or}} & (51) \\{\delta_{app} = {{{1/2}\mspace{11mu}\arcsin\frac{2\; p}{a}} = {{1/2}\mspace{11mu}\arcsin{\frac{2\;{l\left( {v_{l} - v_{r}} \right)}}{b\left( {v_{l} + v_{r}} \right)}.}}}} & (52)\end{matrix}$With δ_(app) one can now extract from (48) the approximate solution

$\begin{matrix}{v_{app} = {\frac{v_{l} + v_{r}}{2}\cos\;\delta_{app}}} & (53)\end{matrix}$and, with (22), from (49)

$\begin{matrix}{h_{app} = {\frac{v_{l} - v_{r}}{2\; b\;\cos\;\delta_{app}}.}} & (54)\end{matrix}$Using (22) one can rewrite equations (47a), (47b) as

$\begin{matrix}{v_{l} = {\frac{v_{app}}{\cos\;\delta_{app}} + {{bh}_{app}\cos\;\delta_{app}}}} & \left( {55a} \right) \\{v_{r} = {\frac{v_{app}}{\cos\;\delta_{app}} - {{bh}_{app}\cos\;\delta_{app}}}} & \left( {55b} \right)\end{matrix}$whence (53) and (54) can also be easily derived.

This approximate solution is usually good enough as long as the wheelangles are not too large, i.e., δ_(app)<δ_(thr), where δ_(thr) is athreshold wheel angle which depends on the geometry of the vehicle andthe precision requirements and will usually be between 30° and 40°.Where δ_(app) is larger than δ_(thr), in particular, if δ_(app)approaches 45°, the approximate solution tends to deviate significantly,yielding often unsatisfactory results. However, it has been found thatthe approximation can be improved to an extent which makes it adequatefor most practical purposes by replacing δ_(app) by δ_(ex) for largevalues of δ_(app) in (53) and (54). This yields results intermediatebetween the exact results according to (41), (44) and the approximateresults according to (53), (54) above which are, however quite close tothe former even for large wheel angles:

$\begin{matrix}{v_{im} = {\frac{v_{l} + v_{r}}{2}\cos\;\delta_{ex}}} & (56) \\{h_{im} = {\frac{v_{l} - v_{r}}{2\; b\;\cos\;\delta_{ex}}.}} & (57)\end{matrix}$

The approximate and intermediate results for the velocity v and theheading rate h differ only by the choice made for a calculation wheelangle δ_(cal) used in the calculations which in the first case equalsδ_(app) and in the second δ_(ex).

δ_(ex) can be calculated directly from δ_(app) as follows:

With

$\begin{matrix}{P = {\frac{v_{l}^{2} - v_{r}^{2}}{v_{l}^{2} + v_{r}^{2}} = \frac{\left( {v_{l} + v_{r}} \right)\left( {v_{l} - v_{r}} \right)}{\frac{1}{2}\left\lbrack {\left( {v_{l} + v_{r}} \right)^{2} + \left( {v_{l} - v_{r}} \right)^{2}} \right\rbrack}}} & (58)\end{matrix}$one gets

${\frac{1}{P} = {{{1/2}\left( {\frac{v_{l} + v_{r}}{v_{l} - v_{r}} + \frac{v_{l} - v_{r}}{v_{l} + v_{r}}} \right)} = {{{1/2}\left( {\frac{1}{p} + p} \right)} = \frac{1 + p^{2}}{2\; p}}}},$i.e.,

$\begin{matrix}{P = \frac{2\; p}{1 + p^{2}}} & (59)\end{matrix}$which, together with (42), leads to

$\begin{matrix}{{{\tan\;\delta_{ex}} = \frac{\frac{2\; p}{a}}{1 + p^{2} + \sqrt{\left( {1 - p^{2}} \right)^{2} - \left( \frac{2\; p}{a} \right)^{2}}}}{and}} & (60) \\{\delta_{ex} = {{\arctan\left( \frac{\frac{2\; p}{a}}{1 + p^{2} + \sqrt{\left( {1 - p^{2}} \right)^{2} - \left( \frac{2\; p}{a} \right)^{2}}} \right)}.}} & (61)\end{matrix}$

Taken together with (51), (61) yields δ_(ex) as a function of δ_(app).

The reverse, i.e., calculation of δ_(app) as a function of δ_(ex), isequally possible. From (20a), (20b) follows

$\begin{matrix}{p = \frac{\begin{matrix}{\sqrt{1 + {2\; a\;\tan\;\delta_{ex}} + {\left( {1 + a^{2}} \right)\tan^{2}\delta_{ex}}} -} \\\sqrt{1 - {2\; a\;\tan\;\delta_{ex}} + {\left( {1 + a^{2}} \right)\tan^{2}\delta_{ex}}}\end{matrix}}{\begin{matrix}{\sqrt{1 + {2\; a\;\tan\;\delta_{ex}} + {\left( {1 + a^{2}} \right)\tan^{2}\delta_{ex}}} +} \\\sqrt{1 - {2\; a\;\tan\;\delta_{ex}} + {\left( {1 + a^{2}} \right)\tan^{2}\delta_{ex}}}\end{matrix}}} & (62)\end{matrix}$and this can be used, together with (52), to calculate δ_(app) fromδ_(ex). In any case, only the fixed parameter a which depends on thegeometry of the vehicle enters into the calculations.

The method according to the invention can be best understood consideringFIG. 2 which provides an overview. For tracking the position of thevehicle, consecutive filter cycles are carried out in the trackingdevice 4, each covering a cycle interval. A complete filter cycle of thetightly coupled Kalman filter which is preferably used comprises aprediction step and, where GNSS measurements are available, a correctionstep which also takes those measurements into account. The use of aKalman filter makes it possible to deal adequately with the randominfluences which usually affect measurements.

A prediction step consists of N time update steps each covering ameasuring interval Δt of, e.g., 0.1 s, with N being, e.g., 10, andconsequently covers a cycle interval NΔt which, in the example, has aduration of 1 s. At the beginning of a prediction step the state vectorof the vehicle is set to the values determined at the end of theprevious cycle interval, or, at the beginning of the process, to someinitial values, which provides a starting point comprising startingvalues x _(ref,0), H₀ for the position and the heading as well as valuesfor the wheel tick calibration factors f_(l), f_(r).

During a time update step the following actions are carried out (s. FIG.3) in the tracking device 4:

A raw velocity T_(l) of the left front wheel 3 a and a raw velocityT_(r) of the right front wheel 3 b are derived from the wheel ticknumbers as registered by the wheel tick units 5 a;b for each of the Ntime update intervals, the numbers being low pass filtered for smoothingout random variations.

Wheel velocities v_(l), v_(r) are calculated from the raw wheelvelocities T_(l), T_(r) assigned to the first measurement interval andfrom the wheel tick calibration factors f_(l), f_(r) which are assumedto be constant throughout the cycle interval using equations (4a), (4b),and new values for the position x _(ref,1) and heading H₁ are determinedwhere equations (1)-(3) with k=0 are used with the calculation wheelangle δ_(cal) having been determined first from the wheel velocitiesv_(l), v_(r) and used to calculate a velocity v and a heading rate h.

k is then (FIG. 2) replaced by k+1 and the process repeated with the endposition from the previous measurement interval used as a starting pointin each case until k=N when the process is terminated after thecalculation of a new state. In this way the trajectory covered by thevehicle during the cycle interval is approximated by a sequence of Nlinear sections each corresponding to a time update step. Thecalculation wheel angle δ_(cal) is assumed to be a constant at everytime update step but the covariance matrix is modified to add processnoise accounting for the random factors which in fact influence thevalue of the said angle used in the calculation.

The result of the prediction step, i.e., of the N^(th) update step,consisting of state and accuracy information, is, where GNSSmeasurements are available, then used in a correction step of thetightly coupled Kalman filter which also processes the said GNSSmeasurements. The filter cycle is finalised and a new state vectordetermined which comprises a corrected position x _(ref,N) and headingH_(N) to be used as a starting point x _(ref,0), H₀ in the first timeupdate step of the prediction step in the next filter cycle as well asrecalibrated wheel tick calibration factors f_(l), f_(r). The positionand heading are output by the tracking device 4 for display or otherpurposes, e.g., for use in processing in a navigation device.

Use of equations (1)-(3) above requires that for every time updateinterval values for v and h are calculated previously which in turnrequires determining the calculation wheel angle δ_(cal) such that v andh can then be calculated as

$\begin{matrix}{v_{im} = {\frac{v_{l} + v_{r}}{2}\cos\;\delta_{cal}}} & (63) \\{h_{im} = \frac{v_{l} - v_{r}}{2\; b\;\cos\;\delta_{cal}}} & (64)\end{matrix}$within a wheel angle interval to which the wheel angle δ is limited andwhich is usually symmetric about 0°, i.e., equals [−δ_(max), +δ_(max)]where δ_(max) is some constant depending on mechanical properties of thevehicle and can, e.g., be somewhere between 35° and 42°. In any case itis assumed that δ_(max) is not greater than δ_(lim).

For the purpose of calculating v and h equation (52) is used in eachcase to determine δ_(app) first. Where δ_(app) is not too large, e.g.,|δ_(app)|≦δ_(thr) where δ_(thr) is a threshold fixed in such a way thatfor |δ_(app)|≦δ_(thr) δ_(app) is within a tolerance range surroundingδ_(ex), δ_(cal) is set equal to δ_(app). The tolerance range isdetermined by the consideration that any angle within the same is closeenough to δ_(ex) to insure that the error introduced by using it instead of δ_(ex) is not significant. In the portion of the wheel angleinterval where δ_(app) is not within the said tolerance range, in theexample in the area defined by |δ_(app)|>δ_(thr), an approximation toδ_(ex) is used which deviates less from the same than δ_(app),preferably a value for δ_(cal) which is virtually equal to δ_(ex). Thatis, δ_(cal) can be determined according to

$\begin{matrix}{\delta_{cal} = \begin{matrix}\delta_{app} & {for} & {{\delta_{app}} \leq \delta_{thr}} \\\delta_{ex} & {for} & {{\delta_{app}} > {\delta_{thr}.}}\end{matrix}} & (65)\end{matrix}$

The value of δ_(cal) from (65) can then be used in (63), (64). Thetolerance range about δ_(ex) can be defined, for instance, by a certainamount or percentage of acceptable deviation, e.g., deviation by notmore than 2°, preferably by not more than 1°. Other criteria, e.g., amaximum acceptable deviation of the velocity v or the heading rate h orboth, in absolute or relative terms, are also possible.

FIG. 4 shows the situation for a=0.31 with an acceptable deviation of±1°. Where δ_(app) is within the tolerance range of δ_(ex)±1° δ_(cal) isset equal to δ_(app). At δ_(app)=δ_(thr) where δ_(app) reaches theboundary of the tolerance range, δ_(cal) is switched to δ_(ex).

The most convenient way to replace δ_(app) by a closer approximation toδ_(ex) if |δ_(app)|>δ_(thr) is to use a look-up table which, for eachone of a row of appropriately spaced values of δ_(ex), contains acorresponding value for δ_(app) such that the values of δ_(app) cover aninterval bounded below by δ_(thr) and above by the value of δ_(app)which corresponds to δ_(ex) when the latter is equal to δ_(max). Anappropriate table can be calculated once for all for a given vehicleusing formulae (52) and (62) and permanently stored. For any value ofδ_(app) with |δ_(app)|>δ_(thr) the corresponding value of δ_(ex) canthen be determined by, e.g., linear, interpolation between neighbouringvalues read from the table. This method does not require much processingpower in the tracking device 4 and memory requirements are also modest.

It is also possible, however, to use a function of the δ's whose valuesdiffer from δ itself, i.e., a function which differs from the identicalfunction.

E.g., from (51) an expression for tan δ_(app) can be derived. Expressingp as

$\begin{matrix}{p = \frac{a\;\tan\;\delta_{app}}{1 + {\tan^{2}\delta_{app}}}} & (66)\end{matrix}$leads to a quadratic equation for tan δ_(app) with the solution

$\begin{matrix}{{\tan\;\delta_{app}} = \frac{1 - \sqrt{1 - \left( \frac{2\; p}{a} \right)^{2}}}{\frac{2\; p}{a}}} & (67)\end{matrix}$which yields, together with (62), tan δ_(app) as a function of tanδ_(ex). In a neighbourhood of zero where this function of p isindefinite it can be replaced by the linear approximation

$\begin{matrix}{{\tan\;\delta_{app}} = {\frac{p}{a}.}} & (68)\end{matrix}$

Again, tan δ_(ex) can vice versa be expressed as a function of tanδ_(app) inserting p from (66) into (60). These relations can be used tocalculate and permanently store a look-up table containing values fortan δ_(ex) as a function of tan δ_(app).

A tan δ_(cal) can then be used in the calculations which is defined as

$\begin{matrix}{{\tan\;\delta_{cal}} = \begin{matrix}{\tan\;\delta_{app}} & {for} & {{\delta_{app}} \leq \delta_{thr}} \\{\tan\;\delta_{ex}} & {for} & {{\delta_{app}} > \delta_{thr}}\end{matrix}} & (69)\end{matrix}$or equivalently

$\begin{matrix}{{\tan\;\delta_{cal}} = \begin{matrix}{\tan\;\delta_{app}} & {for} & {{{\tan\;\delta_{app}}} \leq {\tan\;\delta_{thr}}} \\{\tan\;\delta_{ex}} & {for} & {{{\tan\;\delta_{app}}} > {\tan\;{\delta_{thr}.}}}\end{matrix}} & (70)\end{matrix}$cos δ_(cal) can then be expressed as a function of tan δ_(cal), as

$\begin{matrix}{{{\cos\;\delta_{cal}} = \frac{1}{\sqrt{1 + {\tan^{2}\delta_{cal}}}}},} & (71)\end{matrix}$for calculating v_(im) and −h_(im) according to (63) and (64),respectively.Where (63) and (64) are used it is, however, more convenient to workdirectly with the cosine. An expression for cos δ_(app) can be derivedfrom (51) as follows:p=a√{square root over (1−cos²δ_(app))}cos δ_(app)  (72)leads to a quadratic equation for cos²δ_(app)

$\begin{matrix}{{{\cos^{4}\delta_{app}} - {\cos^{2}\delta_{app}} + \left( \frac{p}{a} \right)^{2}} = 0} & (73)\end{matrix}$with the solution

$\begin{matrix}{{\cos\;\delta_{app}} = {\sqrt{\frac{1 + \sqrt{1 - \left( \frac{2\; p}{a} \right)^{2}}}{2}}.}} & (74)\end{matrix}$for the cosine. Here again, cos δ_(ex) can be expressed as a function ofcos δ_(app), e.g., via the tangent with

${\tan\;\delta_{app}} = {{\frac{\sqrt{1 - {\cos^{2}\delta_{app}}}}{\cos\;\delta_{app}}\mspace{14mu}{and}\mspace{14mu}\cos\;\delta_{ex}} = \frac{1}{\sqrt{1 + {\tan^{2}\delta_{ex}}}}}$and vice versa. This can be used to prepare a look-up table withcorresponding pairs of cos δ_(app) and cos δ_(ex) from which cos δ_(cal)can be easily determined for |δ_(app)|>δ_(thr) or, equivalently, |cosδ_(app)|<cos δ_(thr).

A somewhat different v_(im) can be calculated from tan δ_(ex) using

$\begin{matrix}{v_{im} = {{{{\pm 1}/2}\sqrt{\frac{V_{+}\frac{P}{a}}{\tan\;\delta_{ex}}}} = {{\pm \frac{1}{\sqrt{2\; a}}}\sqrt{\frac{V_{+}p}{\left( {1 + p^{2}} \right)\tan\;\delta_{ex}}}}}} & (75)\end{matrix}$or equivalently

$\begin{matrix}\begin{matrix}{v_{im} = {{\pm \frac{1}{2\sqrt{a}}}\sqrt{\frac{V_{-}}{\tan\;\delta_{ex}}}}} \\{= {{\pm \frac{1}{2\sqrt{a}}}\sqrt{\frac{v_{l}^{2} - v_{r}^{2}}{\tan\;\delta_{ex}}}}} \\{= {{\pm \frac{1}{2\sqrt{a}}}\sqrt{\frac{\left( {v_{l} + v_{r}} \right)\left( {v_{l} - v_{r}} \right)}{\tan\;\delta_{ex}}}}}\end{matrix} & (76)\end{matrix}$which can be easily derived from (41), (42) and yields the exactsolution of equations (20a), (20b) for v. (75) or (76) can, inparticular, be employed where tan δ_(ex) is rather large, e.g., whereδ_(app)>δ_(thr), whereas (63) is used for smaller δ_(app). But it isalso possible to use the said formulae over the whole range, i.e., withtan δ_(cal).h_(im), on the other hand, can also be derived from v_(im)—whethercalculated according to (63) or according to (75) or (76)—and tan δ_(ex)using (22), i.e.,

$\begin{matrix}{{h_{im} = \frac{v_{im}\tan\;\delta_{ex}}{l}},} & (77)\end{matrix}$which, if v_(im) is calculated using (75) or (76), corresponds to

$\begin{matrix}{h_{im} = {{{\pm \frac{1}{2\; l}}\sqrt{V_{+}\frac{P}{a}\tan\;\delta_{ex}}} = {{\pm \frac{1}{\sqrt{2\;{bl}}}}\sqrt{V_{+}\frac{P}{1 + p^{2}}\tan\;\delta_{ex}}}}} & (78)\end{matrix}$or, equivalently, to

$\begin{matrix}{h_{im} = {{{\pm \frac{1}{2\sqrt{bl}}}\sqrt{V_{-}\tan\;\delta_{ex}}} = {{\pm \frac{1}{2\sqrt{bl}}}\sqrt{\left( {v_{l} + v_{r}} \right)\left( {v_{l} - v_{r}} \right)\tan\;\delta_{ex}}}}} & (79)\end{matrix}$which is equivalent to h as calculated according to (44), i.e., theexact solution resulting from equations (20a), (20b). Here again, (77)can be used only where δ_(app)>δ_(thr) or over the whole range, with tabδ_(cal).

Use of the intermediate values v_(im) and h_(im) as calculated accordingto the simple formulae (63) and (64) has been found to be adequate inmost cases even where the calculation wheel angle δ_(cal) varies over aninterval of [−35°, +35°] or more.

If sufficient processing capacity is available, it is also possible touse formulae (52) and (61) in order to calculate δ_(ex) directly fromδ_(app) or formulae (66) and (60) to calculate tan δ_(ex) directly fromtan δ_(app) in stead of using a look-up table or even to use formulae(41) and (44) which yield, together with (23), (24) and (30), a velocityv and heading rate h which are exact solutions of equations (20a) and(20b), corresponding to solutions where δ_(ex) is used for the wheelangle and (75) or (76) and (77) for v and h, respectively.

In the Kalman filter correction step which follows the prediction stepto finalise the filter cycle, apart from GNSS observation equations,observation equations expressing the raw wheel velocities T_(l), T_(r)as functions of the states v, h, f_(l) and f_(r) are used where δ_(cal)figures as a parameter. E.g., equations

$\begin{matrix}{T_{l} = {\frac{1}{f_{l}}\left( {\frac{v}{\cos\;\delta_{cal}} + {{bh}\;\cos\;\delta_{cal}}} \right)}} & \left( {80a} \right) \\{T_{r} = {\frac{1}{f_{r}}\left( {\frac{v}{\cos\;\delta_{cal}} - {{bh}\;\cos\;\delta_{cal}}} \right)}} & \left( {80b} \right)\end{matrix}$are employed which correspond to (55a), (55b) with v_(l), v_(r)expressed according to (4a), (4b), respectively, and δ_(app) replaced byδ_(cal). T_(l) and T_(r) are treated as observations and f_(l), f_(r), vand h as unknowns.

FIG. 5 is a block diagram illustrating one embodiment of the trackingdevice 4. The tracking device 4 comprises an interface 7, a processor 8,and a memory 9. The interface 7 is coupled to the wheel tick units 5 a,5 b by a plurality of input lines and coupled to the antenna forreceiving GNSS signals. The interface 7 is further coupled to theprocessor 8 and the memory 9 for communicating data and commandstherebetween. The interface 7 may also be coupled via a plurality ofoutput lines to a user display (not shown) or other external device toprovide to a user or other device the result of the processing describedherein. The processor 8 and the memory 9 may also be coupled together.The processor 8 performs the processing, calculating, deriving and thelike described above for the tracking device 4 and equations (1) through(80a, 80b).

The memory 9 may function as a data storage system and as a cache memoryfor processing. The memory 9 includes a look-up table 10 that stores thelook-up tables described above.

Some portions of the detailed description above are presented in termsof algorithms and symbolic representations of operations on data bitswithin a computer memory. These algorithmic descriptions andrepresentations are the means used by those skilled in the dataprocessing arts to most effectively convey the substance of their workto others skilled in the art. An algorithm is here, and generally,conceived to be a self-consistent sequence of steps (instructions)leading to a desired result. The steps are those requiring physicalmanipulations of physical quantities. Usually, though not necessarily,these quantities take the form of electrical, magnetic or opticalsignals capable of being stored, transferred, combined, compared andotherwise manipulated. It is convenient at times, principally forreasons of common usage, to refer to these signals as bits, values,elements, symbols, characters, terms, numbers, or the like. Furthermore,it is also convenient at times, to refer to certain arrangements ofsteps requiring physical manipulations of physical quantities as modulesor code devices, without loss of generality.

However, all of these and similar terms are to be associated with theappropriate physical quantities and are merely convenient labels appliedto these quantities. Unless specifically stated otherwise as apparentfrom the above discussion, it is appreciated that throughout thedescription, discussions utilizing terms such as “processing” or“computing” or “calculating” or “determining” or “displaying” or thelike, refer to the action and processes of a computer system, or similarelectronic computing device, that manipulates and transforms datarepresented as physical (electronic) quantities within the computersystem memories or registers or other such information storage,transmission or display devices.

Certain aspects of the present invention include process steps andinstructions described herein in the form of an algorithm. It should benoted that the process steps and instructions of the present inventioncould be embodied in software, firmware or hardware, and when embodiedin software, could be downloaded to reside on and be operated fromdifferent platforms used by a variety of operating systems.

The present invention also relates to an apparatus for performing theoperations herein. This apparatus may be specially constructed for therequired purposes, or it may comprise a general-purpose computerselectively activated or reconfigured by a computer program stored inthe computer. Such a computer program may be stored in acomputer-readable storage medium, such as, but is not limited to, anytype of disk including floppy disks, optical disks, CD-ROMs,magnetic-optical disks, read-only memories (ROMs), random accessmemories (RAMs), EPROMs, EEPROMs, FLASH, magnetic or optical cards,application specific integrated circuits (ASICs), or any type of mediasuitable for storing electronic instructions, and each coupled to acomputer system bus. Furthermore, the computers referred to in thespecification may include a single processor or may be architecturesemploying multiple processor designs for increased computing capability.

The algorithms and displays presented herein are not inherently relatedto any particular computer or other apparatus. Various general-purposesystems may also be used with programs in accordance with the teachingsherein, or it may prove convenient to construct more specializedapparatus to perform the required method steps. The required structurefor a variety of these systems will appear from the description above.In addition, the present invention is not described with reference toany particular programming language. It will be appreciated that avariety of programming languages may be used to implement the teachingsof the present invention as described herein, and any references belowto specific languages are provided for disclosure of enablement and bestmode of the present invention.

In addition, the language used in the specification has been principallyselected for readability and instructional purposes, and may not havebeen selected to delineate or circumscribe the inventive subject matter.Accordingly, the disclosure of the present invention is intended to beillustrative, but not limiting, of the scope of the invention, which isset forth in the claims.

What is claimed is:
 1. A method for tracking a position and a heading ofa vehicle, at least intermittently using dead reckoning tracking basedon velocities of left and right front wheels, the method comprising:receiving, by using a processor, measurement parameters indicative of avelocity v_(l) of the left front wheel and measurement parametersindicative of a velocity v_(r) of the right front wheel, deriving, byusing the processor, a value for the velocity v_(l) of the left frontwheel from the parameters indicative of a velocity v_(l) of the leftfront wheel and a value of the velocity v_(r) of the right front wheelfrom the parameters indicative of a velocity v_(r) of the right frontwheel, calculating, by using the processor, a function of a calculationwheel angle δ_(cal) from the velocity v_(l) of the left front wheel andthe velocity v_(r) of the right front wheel, the calculation wheel anglevarying over a wheel angle interval, using, by the processor, thefunction of the calculation wheel angle δ_(cal) to calculate a velocityv and a heading rate h of the vehicle, where the calculation wheel angleδ_(cal) reflected in said function deviates at least in a portion of thewheel angle interval, from an exact wheel angle δ_(ex) calculatedaccording to$\delta_{ex} = {\arctan\left( \frac{\frac{2\; p}{a}}{1 + p^{2} + \sqrt{\left( {1 - p^{2}} \right)^{2} - \left( \frac{2\; p}{a} \right)^{2}}} \right)}$by less than an approximate wheel angle defined as$\delta_{app} = {{1/2}\mspace{11mu}\arcsin\frac{2\; p}{a}}$ deviatesfrom the same, where a is one half of the front track gauge divided bythe wheel base and ${p = \frac{v_{l} - v_{r}}{v_{l} + v_{r}}},$calculating, by using the processor, values of the position and theheading of the vehicle based on a previous position and heading of thevehicle, using the calculated velocity v and heading rate h of thevehicle.
 2. The method according to claim 1, wherein for the whole wheelangle interval, the extent of the deviation of the calculation wheelangle δ_(cal) from the exact wheel angle δ_(ex) is smaller than 2°. 3.The method according to claim 1, further comprising calculating afunction of an approximate wheel angle δ_(app), and, if the approximatewheel angle δ_(app) is smaller than a threshold wheel angle, thefunction of the wheel angle used in the calculation of the heading rateh reflects that the calculation wheel angle δ_(cal) equals theapproximate wheel angle δ_(app) whereas, if the approximate wheel angleδ_(app) is greater than the threshold wheel angle, the function reflectsthat the calculation wheel angle δ_(cal) differs from the approximatewheel angle δ_(app) and is determined from the function of the latter.4. The method according to claim 3, where for the determination of thefunction of the calculation wheel angle δ_(cal) from the function of theapproximate wheel angle −δ_(app) a look-up table is used where pairsconsisting of the value of the function of the exact wheel angle δ_(ex)and the function of the corresponding value of the approximate wheelangle δ_(app) are stored.
 5. The method according to claim 4, furthercomprising calculating, using interpolation, the value of the functionof the exact wheel angle δ_(ex), where the value of the function of theapproximate wheel angle δ_(app) falls between two values of the samestored in the look-up table.
 6. The method according to claim 1, whereinthe velocity v of the vehicle is calculated as$v = {\frac{v_{l} + v_{r}}{2}\cos\;{\delta_{cal}.}}$
 7. The methodaccording to claim 1, wherein the velocity v of the vehicle iscalculated as$v = {\frac{1}{2\sqrt{a}}{\sqrt{\frac{v_{l}^{2} - v_{r}^{2}}{\tan\;\delta_{cal}}}.}}$8. The method according to claim 1, wherein the heading rate h iscalculated as $h = \frac{v\;\tan\;\delta_{cal}}{l}$ where l is the wheelbase.
 9. The method according to claim 1, wherein the heading rate h ofthe vehicle is calculated as$h = \frac{v_{l} - v_{r}}{2\; b\;\cos\;\delta_{cal}}$ where b is onehalf of the front track gauge.
 10. The method according to claim 1,wherein the calculation wheel angle δ_(cal) varies over a rangecomprising the interval [−35°, 35° ].
 11. The method according to claim1, further comprising filtering in consecutive filter cycles, whereineach filter cycle comprises a prediction step consisting of at least onetime update step where in each case the values of the position and ofthe heading of the vehicle are calculated from the velocity v and theheading rate h, each calculation being based on a previously establishedvalue of the position and the heading, and, if Global NavigationSatellite System measurements are available, a correction step where thevalues of the position and the heading calculated in the prediction stepare corrected with the Global Navigation Satellite System measurementsbeing taken into account.
 12. The method according to claim 11, whereinevery prediction step comprises a fixed number of consecutive timeupdate steps.
 13. The method according to claim 11, wherein at eachcorrection step wheel calibration factors which are used to determinethe velocity v_(l) of the left front wheel and the velocity v_(r) of theright front wheel from raw wheel velocities which are based on measuredwheel tick numbers are recalibrated.
 14. A tracking device usable in avehicle, comprising: a plurality of input lines for receiving data fromwheel tick units; a plurality of output lines for providing at least areference position; a data storage system; a processor coupled to theplurality of input lines, the plurality of output lines, and the datastorage system and configured to receive sensor measurements indicativeof the velocity v_(l) of a left front wheel and sensor measurementsindicative of the velocity v_(r) of a right front wheel, derive a valuefor the velocity v_(l) of the left front wheel from the sensormeasurements indicative of a velocity v_(l) of the left front wheel anda value of the velocity v_(r) of the right front wheel from the sensormeasurements indicative of a velocity v_(r) of the right front wheel,calculate a function of a calculation wheel angle δ_(cal) from thevelocity v_(l) of the left front wheel and the velocity v_(r) of theright front wheel, the calculation wheel angle δ_(cal) varying over awheel angle interval, use the function of the calculation wheel angleδ_(cal) to calculate a velocity v and a heading rate h of the vehicle,where the calculation wheel angle δ_(cal) which is reflected in saiddeviates, at least in a portion of the wheel angle interval, from anexact wheel angle δ_(ex) calculated according to$\delta_{ex} = {\arctan\left( \frac{\frac{2\; p}{a}}{1 + p^{2} + \sqrt{\left( {1 - p^{2}} \right)^{2} - \left( \frac{2\; p}{a} \right)^{2}}} \right)}$by less than an approximate wheel angle defined as$\delta_{app} = {{1/2}\mspace{11mu}\arcsin\frac{2\; p}{a}}$ deviatesfrom the same, where a is one half of the front track gauge divided bythe wheel base and ${p = \frac{v_{l} - v_{r}}{v_{l} + v_{r}}},$calculate values of a position and a heading of the vehicle based on aprevious position and heading of the vehicle, using the calculatedvelocity v and heading rate h of the vehicle.
 15. The tracking deviceaccording to claim 14, wherein for the whole wheel angle interval, theextent of the deviation of the calculation wheel angle δ_(cal) from theexact wheel angle δ_(ex) is smaller than 2°.
 16. The tracking deviceaccording to claim 14, wherein the processor is further configured tocalculate a function of the approximate wheel angle δ_(app) and, if theapproximate wheel angle δ_(app) is smaller than a threshold wheel angle,to use the function of the wheel angle in the calculation of the headingrate h which reflects that the calculation wheel angle δ_(cal) equalsthe approximate wheel angle δ_(app) whereas, if the approximate wheelangle δ_(app) is greater than the threshold wheel angle, the functionreflects that the calculation wheel angle δ_(cal) differs from theapproximate wheel angle δ_(app) and is determined from the function ofthe latter.
 17. The tracking device according to claim 16, wherein theprocessor is further configured to use, for the determination of thefunction of the calculation wheel angle δ_(cal) from the function of theapproximate wheel angle δ_(app), a look-up table where pairs consistingof the value of the function of the exact wheel angle δ_(ex) and thefunction of the corresponding value of the approximate wheel angleδ_(app) are stored.
 18. The tracking device according to claim 17,wherein the processor is further configured to calculate, usinginterpolation, the value of the function of the exact wheel angleδ_(ex), where the value of the function of the approximate wheel angleδ_(app) falls between two values of the same stored in the look-uptable.
 19. The tracking device according to claim 14, having at leastone further input line for receiving a signal from a Global NavigationSatellite System antenna, the processor being configured to carry outconsecutive filter cycles, wherein each filter cycle comprises aprediction step consisting of at least one time update step where ineach case the values of the position and the heading of the vehicle arecalculated from the velocity v and the heading rate h, each calculationbeing based on a previously established value of the position and theheading, and, if Global Navigation Satellite System measurements can bederived from the signal received from the Global Navigation SatelliteSystem antenna, a correction step where the values of the position andthe heading calculated in the prediction step are corrected with theGlobal Navigation Satellite System measurements being taken intoaccount.